﻿#include "relaxedJacobiSmoother.hpp"

template <typename Cmpt>
void UNAP::RelaxedJacobiSmoother<Cmpt>::smooth(Vector<Cmpt> &x,
                                               const Matrix<Cmpt> &A,
                                               const Vector<Cmpt> &b,
                                               const label nSweeps) const {
  const label rowSize = x.size();
  Cmpt *xPtr = x.data();
  const Cmpt *diagPtr = A.diag().data();
  const Cmpt *bP = b.data();
  Vector<Cmpt> tmp(rowSize);
  Cmpt *tmpP = tmp.data();

  for (label sweep = 0; sweep < nSweeps; sweep++) {
    A.spMV(tmp, x);
    forAll(cellI, rowSize) {
      tmpP[cellI] = -(omega_ / diagPtr[cellI]) * (tmpP[cellI] - bP[cellI]);
      xPtr[cellI] += tmpP[cellI];
    }
  }
}

template class UNAP::RelaxedJacobiSmoother<scalar>;
